Lab5

Task 1

#create GEOID column in bmps from subset GeographyName
bmps$GEOID = stringr::str_sub(bmps$GeographyName, 1, 5)

#sum cost of bmps per geoid
tot.cost <- bmps %>% filter(., Cost > 0) %>% #filter out values of 0
  group_by(GEOID) %>% #GEOID to get a single value to join to counties later
  summarise(., TotalBMPCosts = as.numeric(format(round(sum(Cost),2), nsmall=2))) #trims costs to 2 decimal places to match currency

#join tot.cost to counties on GEOID, drop counties with no BMP costs
bmp.sf <- left_join(CBWcounties, tot.cost, by = c("GEOID10"="GEOID")) %>% 
  filter(., TotalBMPCosts > 0) %>% 
  st_transform(., crs = "EPSG:4326")

#create bins and palette
maxB  <- max(bmp.sf$TotalBMPCosts) + 1 #find max value and add 1 for upper bin
bins <- c(0, 0.2*maxB, 0.4*maxB, 0.6*maxB, 0.8*maxB,maxB) #create quantile classes programmaticaly
pal <- colorBin("BuPu", domain = bmp.sf$TotalBMPCosts, bins=bins) #create palette

task1 <- leaflet(bmp.sf) %>% 
  addTiles() %>% 
  addPolygons(fillColor = ~pal(TotalBMPCosts),
              weight = 2,
              opacity = 1,
              color = "blue",
              fillOpacity = 1,
              label = paste0("$", bmp.sf$TotalBMPCosts)) %>% 
  addLegend(pal = pal,
            values = ~TotalBMPCosts,
            title = "Total BMP Costs", 
            position = "bottomright")

task1

Task 2

#percentage of households that are single parent
aoi <- aoi %>% mutate(., oneparent = ((DP0130007 + DP0130009)/ DP0130001) * 100) #single male and female parents
                                                                           # over total households
#spatial weights matrix
aoi.lw <- poly2nb(aoi, queen = TRUE) %>% nb2listw(., style = "W")

#calculate LISA
aoi$LI <- localmoran(aoi$oneparent, aoi.lw)[, 1]
aoi$pval <- localmoran(aoi$oneparent, aoi.lw)[, 5]

#covert aoi to WGS84, otherwise Leaflet will show single worldwide polygon
aoi.wgs <- st_transform(aoi, crs = "EPSG:4326")

#bins and palette
max2 = max(aoi.wgs$LI) + 0.1
min2 = min(aoi.wgs$LI) - 0.1
dff2 = max2 - min2
bins2 <- c(min2, max2-(0.75*dff2), max2-(0.5*dff2), max2-(0.25*dff2), max2)
pal2 <- colorBin("Reds", domain = aoi.wgs$LI, bins=bins2)

#make map
task2 <- leaflet(aoi.wgs) %>% #all data come from aoi.wgs
  addProviderTiles(providers$Esri.WorldGrayCanvas, group = "Gray Canvas (Default)") %>% #default gray canvas
  addProviderTiles(providers$Stamen.TerrainBackground, group = "Stamen Terrain") %>%  #terrain
  addProviderTiles(providers$OpenStreetMap, group = "OSM") %>% #street map if needed
  addPolygons(fillColor = ~pal2(LI),
              weight = 1,
              opacity = 1,
              color = "white",
              fillOpacity = 0.9, 
              label = paste0("p-value: ", aoi$pval)) %>% 
  addLegend(pal = pal2,
            values = ~LI,
            title = "LISA", 
            position = "bottomright") %>% 
  addLayersControl(
    baseGroups = c("Gray Canvas (Default)", "Stamen Terrain", "OSM"), #assign basegroups to provider tiles above
    options = layersControlOptions(collapsed = FALSE)
  )
  
task2

#bins and palette
max3 = max(lnk.tracts$pov.perc) + 0.1 #define max
min3 = min(lnk.tracts$pov.perc) - 0.1 #define min
dff3 = max3 - min3 #find difference
bins3 <- c(min3, max3-(0.75*dff3), max3-(0.5*dff3), max3-(0.25*dff3), max3) #4 classes using difference ratios
pal3 <- colorBin("YlOrRd", domain = lnk.tracts$pov.perc, bins=bins3)

#steams wgs84
streams.wgs <- streams %>% st_simplify(.) %>% #needed to simplify for leaflet
  st_transform(., crs = "EPSG:4326")

streams.wgs <- streams.wgs[lnk,]

# make maps
task3 <- leaflet() %>% #mapping more than one important data source so not added
  addProviderTiles(providers$Esri.WorldGrayCanvas, group = "Gray Canvas (Default)") %>% 
  addProviderTiles(providers$OpenStreetMap, group = "OSM") %>%
  addRasterImage(dem, #add dem
                 group = "DEM") %>% 
  addPolygons(data = lnk.tracts, #add tracts and bins created above
              fillColor = ~pal3(pov.perc),
              weight = 1,
              opacity = 1,
              color = "gray",
              fillOpacity = 0.9,
              group = "Tracts") %>% 
  addPolygons(data = lnk.parks, #add parks with green color
              fillColor = "#098908",
              weight = 0.25,
              opacity = 1,
              color = "black",
              fillOpacity = 0.9,
              group = "Parks") %>% 
  addPolygons(data = lnk.lakes, #add lakes 
              fillColor = "blue",
              weight = 0.25,
              opacity = 1,
              color = "black",
              fillOpacity = 0.9,
              group = "Water") %>% #water group for lakes and streams
  addPolylines(data = streams.wgs, #add streams
               color = "blue",
               weight = 2, #weight of 2 is good blend of visibility and not overdoing it
               opacity = 1,
               group = "Water") %>% #water group for lakes and streams
  addLegend(data = lnk.tracts, #add legend for data of interest
            pal = pal3,
            values = ~pov.perc,
            title = "% in Poverty",
            position = "bottomright") %>%
  addLayersControl(
    baseGroups = c("Gray Canvas (Default)" ,"OSM"),
    overlayGroups = c("DEM", "Tracts", "Parks", "Water"),
    options = layersControlOptions(collapsed = FALSE)) %>% 
  addMeasure() #add ability to measure distances between 

## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj =
## prefer_proj): Discarded ellps WGS 84 in Proj4 definition: +proj=merc +a=6378137
## +b=6378137 +lat_ts=0 +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null
## +wktext +no_defs +type=crs

## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj =
## prefer_proj): Discarded datum World Geodetic System 1984 in Proj4 definition

## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj =
## prefer_proj): Discarded ellps WGS 84 in Proj4 definition: +proj=merc +a=6378137
## +b=6378137 +lat_ts=0 +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null
## +wktext +no_defs +type=crs

## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj =
## prefer_proj): Discarded datum World Geodetic System 1984 in Proj4 definition

task3

Questions

1. Reflect on the labs from this semester. What did you learn? What did you like? What did you not

like?

I learned an awful lot this semester. While I had previous experience with R, I didn’t truly understand it until I got an opportunity to learn it in the spatial context. In particular, using the sf package and tidyverse together made spatial data manipulation so much faster than many of the GUI based GIS programs. In addition, I felt like I understand spatial analysis better and being able to script all of these steps made it stick. I really enjoyed the web mapping lab because the leaflet library makes mapping a cinch. ALso, cleaning and merging data in the earlier labs made the process really smooth the rest of the semester. I will say the labs did take longer to complete than I thought they would, even though I = had some previous R and programming experience. This could be a barrier for future students with less experience than myself. In addition, I felt the course did not have as much spatial analysis or advanced spatial analysis (autoregressions, etc) as I would have hoped but for learning GIS and basic SA it fit the bill perfectly.

2. Describe the “one thing” you chose to add to your map in Task 3 above. What did you do, and why

is it applicable to your map?

I chose to add the measuring tool so that if a user can measure distances from a respective location to the nearest park in order to understand poverty and the lack of equitable access to green spaces such as parks.